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Abstract. We consider the numerical approximation of a general second order semi— linear 
parabolic stochastic partial differential equation (SPDE) driven by additive space-time noise. We 
introduce a new scheme using in time a linear functional of the noise with a semi— implicit Euler— 
Maruyama method and in space we analyse a finite element method although extension to finite 
differences or finite volumes would be possible. We give the convergence proofs in the root mean 
square norm for a diffusion reaction equation and in root mean square norm in the presence 
of advection. We examine the regularity of the initial data, the regularity of the noise and errors 
from projecting the noise. We present numerical results for a linear reaction diffusion equation in two 
dimensions as well as a nonlinear example of two-dimensional stochastic advection diffusion reaction 
equation. We see from both the analysis and numerics that we have better convergence properties 
over the standard semi-implicit Euler-Maruyama method. 

Key words. Parabolic stochastic partial differential equation, Finite element. Modified Semi- 
implicit Eulcr-Maruyama, Strong numerical approximation. Additive noise. 

1. Introduction. We analyse the strong numerical approximation of Ito stochas-| 

tic partial differential equation defined in 51 C . Boundary conditions on the domain 
ri are typically Neumann, Dirichlet or some mixed conditions. We consider 

dx ^ {AX + F{x))dt + dw, x{Q)^Xo, te[o,r], T>0 (1.1) 

in a Hilbert space H = L'^{VL). Here A is the generator of an analytic semigroup 
S{t) := e*-^ , t > with eigenfunctions and eigenvalues A^, i G N''. F is a nonlinear 
function of X and possibly VX and the noise term W{x,t) is a Q- Wiener process 
that is white in time and defined on a filtered probability space (D, J^, P, {.Fi }j>q). 
The noise can be represented as 

W{x,t)^J2VQ^e^ix)P^it), (1.2) 

where qi > 0, i E and /Si are independent and identically distributed standard 
Brownian motions. Precise assumptions on A, F and W are given in Section [2] and 
under these type of technical assumptions it is well known, see |5l UHl |4] that the 
unique mild solution is given by 

X{t)^S{t)Xo+ f S{t- s)F{X{s))ds + 0{t) (1.3) 
Jo 

with the stochastic process O given by the stochastic convolution 

0{t) = f S{t-s)dW{s). (1.4) 

^Department of Mathematics and Maxwell Institute, Heriot-Watt University, Edinburgh EH14 
4AS, UK 

-fg . j . lordOhw .ac.uk 

§atl500hw. ac .uk (Antoine.Tambue@math.uib.no). Corresponding author. 

'Department of Mathematics, University of Bergen, Johannes Bruns Gate 12, 5020 Bergen .Nor- 
way. 



2 



G.J. Lord and A. Tambue 



Typical examples of the above type of equation are stochastic (advcction) reaction 
diffusion equations where A = A arising from example in pattern formation in physics 
and mathematical biology. We illustrate our work with both a simple reaction diffu- 
sion equation where we can construct an exact solution 

dX = (DAX -\X)dt + dW (1.5) 

as well as the stochastic advection reaction diffusion equation 

dX ^ (^DAX -\J -{ciX)^ j^^^^-^dt + dW (1.6) 

where D > is the diffusion coefficient, q is the Darcy velocity field [5] and A is a 
constant depending of the reaction. 

The study of numerical solutions of SPDEs is an active area of research and there is 
a growing literature on numerical methods for SPDEs. We analyse convergence of a 
finite element discretization in space combined with a semi-implicit discretization in 
time that uses a linear functional to approximate the noise. This approach extends 
the analysis of Jentzen [Til SHI HI] which is based on a Fourier spectral discretization. 
This discretization diagonalizes the linear operator A and then exploits that in each 
Fourier mode the noise is an Ornstein-Uhlenbeck process and hence the variance 
is known. For complex domains or mixed boundary conditions a Fourier spectral 
approach is not feasible and, for example, finite element (finite difference or finite 
volume) discretization is preferred even they are less accurate compared to spectral 
methods, however this destroys the diagonalization of the linear operator. 
Note that even if the eigenfunctions of the linear operator are known explicitly for 
the given domain and boundary conditions it may be that a finite element (or finite 
volume) discretization is prefered over a Fourier spectral method. For example where 
the solution has steep fronts finite elements and finite volumes are efficient and there is 
no issue with any Gibbs phenomena as in Fourier spectral approaches. Such solutions 
arise, for example, in heterogeneous porous media. In this case the Darcy velocity q 



in (1.6) is solved for from the following system 



Vq = 

where k is the permeability tensor, p is the pressure and /i is the dynamic viscosity 
of the fluid ^2 . Typically this would be solved using finite elements or finite volumes. 
Furthermore, due to the heterogenous nature of the permeability such problems often 
naturally give rise to non-uniform grids. Our approach, with a projection of the noise 
onto a standard flnite element grid, allows practitioners to simply adapt their existing 
codes to examine the effects of stochastic forcing. 

We perform our analysis for the case of finite elements and numerically we also exam- 
ine finite volumes. Our work differs from other finite element discretizations [T1 [H1[T5] 
where the approach to the noise is to consider it directly in the finite element space. 
We follow more closely [22l ESI Ej and introduce a projection onto a finite number 
modes and a projection onto the finite element space. The aim is to gain the flexibil- 
ity of the flnite element (finite volume) discretization to deal with flow and transport 



problems (1.61- (1.7), complex boundaries, mixed boundary conditions and inhomo- 



geneous boundary conditions as well as the overcoming the order barrier as in 
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We also show in Section [4~2] that equipped with the eigenvalues and eigenfunctions of 
the operator A with Neumann or Dirichlet boundary conditions, we can apply the new 
scheme with mixed boundary conditions for the operator A = DA without explicitly 
having the eigenvalues and eigenfunctions of A. In general for an operator A with 
given eigenvalues and eigenfunctions for Neumann or Dirichlet boundary conditions, 
the new scheme can be used for the operator A with mixed bounded conditions without 
explicit eigenvalues and eigenfunctions. 

We give convergence proofs in root mean square L^(f2) norm for reaction-diffusion 
equations and in root mean square H^{fl) norm for advection reaction-diffusion for 
spatially regular noise. The smoothing effect of the semigroup generated by the 



operator A in the SPDE (1.1 1 and various semigroup estimates play an important 
role in the proofs. 

The paper is organised as follows. In Section [2] we present the numerical scheme and 
assumptions that we make on the linear operator, nonlinearity and the noise. We 
then state and discuss our main results. In Section |3.2| and Section |3.3| we present 
the proof of the convergence theorems. We end by presenting some simulations in 
Section |4j these are applied both to a linear example where we can compute an exact 
solution as well as a more realistic model coming from model of the advection and 
diffusion of a solute in a porous media with a non-linear reaction term. 

2. Numerical scheme and main results. Let us start by presenting briefly 
the notation for the main function spaces and norms that we use in the paper. We 
denote by || • || the norm associated to the inner product (•, •) of the K— Hilbert space 
H = L^(r2). For a Banach space V we denote by |j • ||v the norm of the space V, L{V) 
the set of bounded linear mapping from V to V, i*^^-'(V) the set of bounded bilinear 
mapping from V x V to M and i2(D, V) the space defined by 

L2(D, V) = |u random variable with value in V : EjlwHy = J ||u(aj)||y(iP(aj) < oo| . 

Throughout the paper we assume that fl is bounded and has a smooth boundary or 
is a convex polygon. For convenience of presentation we take A to be a self adjoint 
second order operator as this simplifies the convergence proof. More precisely 



^ = V.DV(.)+Do,oI= 2^^(A,,^l+I?o,oI, (2.1) 

where we assume that Di j £ L°°{H), and that there exists a positive constant ci > 
such that 



d 

^ A,j(a;)6^j > cilfp V^eM'* xeU ci > 0. (2.2) 



Under condition ( |2.2[ ), it is well known (see [7]) that the linear operator A generates 
an analytic semigroup S{t) = e*"*. We introduce two spaces H and V where H C F, 
that depend on the choice of the boundary conditions for the domain of the operator 
A and for the variational form associated to the operator A. For Dirichlet boundary 
conditions we let 

V = M^ H^in) ^{veH^{n):v = on dfl}. 
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For Robin boundary conditions (Neumann boundary condition being a particular 
case) V = H^{n) and 

m={veH^{n):dv/diyA + crv^O on dn} , a eR. 

Functions in H can satisfy the boundary conditions and with H in hand we can char- 
acterize the domain of the operator (— yl)''/^ and have the following norm equivalence 
[161 [6] for r e {1,2} 

MH^^n) ^ II (-^)''/''^ll =: Ikll. V « e v{{-AY^^) = MnH^n). 

In the Banach space V{{-A)"'^^), a S E, we use the notation := \\{-A)"^^.\\. 
We recall some basic properties of the semi group S{t) generated by A. 
Proposition 2.1. [Smoothing properties of the semi groupf^] 

Let a > 0, P >Q and < 7 < 1, then there exist C > such that 

\\{-AfS{t)\\L(LHm<Ci'" V t>0 
\\{-A)-\I~-Sm\LiLHn))<Ct'^ for t>0. 

In addition, 

{-AfS{t) = S{t)i-A)^ on V{{~Af) 

If P>1 then V{{^A)f^) c V{{-Ay), 
\\D\S{t)v\\p<Ct-'-^^-^y^\\v\U, t > 0, t; e Z = 0,1, 

where D\ := — j . 

We consider discretization of the spatial domain by a finite element triangulation. Let 
7?i be a set of disjoint intervals of $7 (for d = 1), a triangulation of (for d = 2) or a set 
of tetrahedra (for d — i). Let Vh C V denote the space of continuous functions that 
are piecewise linear over the triangulation Th- To discretize in space we introduce two 
projections. Our first projection operator is the L^(f2) projection onto Vh defined 
for u e L^(ri) by 

{PhU,x)^iu,x) Vxen. (2.3) 
Then Ah ■ Vh — >■ Vh is the discrete analogue of A defined by 

{Ahip, x) - {Aif, x) X e Vh. (2.4) 

We denote by Sh the semigroup generated by the operator Ah- 

The second projection Pn, e N is the projection onto a finite number of spectral 
modes defined for u G L^(il) by 

Pnu = ^ {ei,u)ei, 



where Xn = {1,2,..., N} . The semi-discrete in space version of the problem ( 1.1 ) is 
to find the process X''{t) = X''{.,t) e Vh such that for t G [0,T], 

dX'' = (AhX'' + PhF{X''))dt + PhPNdW, X'^iO) = PhXa. (2.5) 
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We denote by X the solution of the random system 



X\t) = Shm^^iO) + / Sh{t - s)F{X''{s))ds 



As in (1.3 1, by splitting we have 

X^{t) = x\t) + PhPNO{t). 

We now discretize in time by a semi-impUcit method to get the fully discrete approx- 
imation of X defined by Z'^ 



7n — 1 



k=0 



(2.6) 



where 

It is straightforward to show that 

ShAt {Zln + ^tPhF{Zt + PhPNOitm))) 



(2.7) 



(2.8) 



Finally we can define our approximation X!^ to X{tm), the solution of equation (1.1) 
by 



PhPNO{U^). 



(2.9) 



Therefore 



ShAt {Xt + At P,,F(X^) - PhPNO{t^)) + PhPNO{Ur^+i), (2.10) 



where according to (1.4 1, we generate 0(<m+i) from 0(i,n) by 



0(Wi) 



AtA 



0% 



-s)A 



dW{s). 



(2.11) 



If we assume that Q has the same eigenfunctions as the linear operator A and we 
have diagonalized the operator A by a Fourier spectral method then (2.111 reduces 
to an Ornstein-Uhlenbeck process in each Fourier mode as in [T^. This can then be 
simulated numerically with the correct mean and variance in each mode. The new 
scheme (2.10) uses a finite element discretization and projects the linear functional 
of the noise and hence we expect superior approximation properties over a standard 
semi-implicit Euler discretization for a finite element discretization. 
The existence and uniqueness of the solution of equation is well known fjl3 | [5l[T8 l l4]\ 
For convergence proofs below we need sufficient regularity of the mild solution X, 
and therefore we will use some weak assumptions. We describe in detail the weak 
assumptions that we make on the linear operator A, the nonlinear term F and the 
noise dW . 

Assumption 2.2. [Linear operator] The linear operator —A given in (2.1) is 
positive definite. Then there exists sequences of positive real eigenvalues {A„}„gpfd 
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with inf Ai > and an orthonormal basis in H of eigenf unctions {eijigNd such that 
the linear operator —A : 2?(— A) <Z H ^ H is represented as 

-Av = ^ A,(e„ v)e, V t> G V{-A) 

where the domain of —A, A) = {v E H : ^ A^|(ei,w)| < oo}. 
Assumption 2.3. [Nonlinearity] 

We assume that there exists a positive constant L > such that F satisfies one of the 
following. 

(a) Let V be a separable Banach space such that 2?((-A)i/2) c V C = L'^{n) 
continuously. 

-F ; V — !■ V is a twice continuously Nemytskii Frechet differentiable mapping with 
\\F'iv)w\\ < L\\w\\, ||F^'(«)|U(V) < L, < L 

and 

<i(i + ll"llp((-A)i/=)) V v,wev, ueVii^A)'/^), 

where {F'{u))* is the adjoint of F'{u) defined by 

{{F'{u))*v,w) = {v,F'{u)w) V v,weH^L^{Q). 
As a consequence 

<L||z-y|| yz,YeH, 

andyy eH = L'^in) 

\\F{Y)\\ < ||F(0)|| + \\F{Y) - ^^(0)11 < ||F(0)|| + L\\Z\\ < C(||F(0)|| + 

(b) F is globally Lipschitz continuous from (_ff-'^(f2), ||.||^i(q)) to {H — i^(f2), ||.||), 



I.e. 



\\F{Z) - F{Y)\\ <L\\Z- VZ, Y e H\n). 



Remark 2.4. As the M— Hilbert space is H = L'^{Vl), F'iu) is self-adjoint we have 
F'{u)* — F'[u), u G !?((— A)"'^/^). This come obviously to the fact that 



{v,w) 



/ vwdx v,w & L^{rt). 



We keep the adjoint notation F' [u)* in Assumption 2.3 as results in this paper can 



be generalized in any Hilbert space H. Possibly expression of F is 

F{u){x) ^ f{x,u{x)) xeQ, 

where / : f7 x M — > M is twice continuously differentiable function with the bounded 
partial derivatives. 
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We note that in condition (a) of Assumption 2.3. the possible choice of V can be H 
or H^{ fl). Strong convergence proofs with a Nemytskii operator for the nonlinearity 
were recently obtained in 

We now turn our attention to the noise term and introduce spaces and notation that 
we need to define the Q- Wiener process. An operator T G L{H) is Hilbert-Schmidt if 

mis ■■=J2m^\\'<^- 

The sum in || -Hl^g is independent of the choice of the orthonormal basis in H. We de- 
note the space of Hilbert- Schmidt operator from Q^/'^{H) to H by L§ := HS{Q^/'^{H), 
and the corresponding norm 1 1 . 1 1 by 

\\T\\,.:^\\TQ''^\\ns-iY.\\TQ"'^A T ^ L\. 

Let be a continuous L2~P^ocess. We have the following equality known as the Ito's 
isometry 

E|| / ipdWf= f EMlods^ f E||(^gi/2||2^5ds. 
Jo Jo ^ Jo 

We assume sufficient regularity of the noise for the existence of a mild solution and to 
project the noise into the space Vh ■ To be specific we need that the stochastic process 
O is in or in space. This is used in our proofs to apply the finite clement 
projection P/j in the errors estimates. 

Notice that for all t g [0,T] the process 0{t) is an adapted stochastic process to the 
filtration {J-t)t>o with continuous sample paths such that 0{t2) — S{t2~ti)0{ti), < 
ti < t2 < T is independent of J^t^ . 

Assumption 2.5. [Noise and stochastic process O] We assume one or both of 
the following. 

(a) The sequence [qi) in the noise representation (|1.2[) satisfies 



J2 Ap'ft <oo^\\ (-A)('-i)/2gi/2||^^ < oo r e {1, 2} . 
As a consequence 

0{t) e L2{D,V{{^AY^^)), 0<t<T (2.12) 

v{{-Ay/^) = H n H'-in), r e {i, 2} . (2.13) 

(b) For some 6 £ (0, 1/2] and positive constant C > 

E (||0(<2) - O(ii)llt) < C{t2 - t^f' 0<ti<t2<T, 



where V is defined in Assumption 2.3(a). 
We note again that we have assumed only for convenience that Q and A have the 
same eigenfunctions. In particular it is useful to characterize the decay of qi in the 
condition (a) of Assumption 2.3 above. 



One can prove that in condition (b) of Assumption 2.5 for V = H = L {Q) we can 
takee = 1/2. For V = H^{n) we can take = l/2if 0(<) e L2(D, 2?(-A)), 0<t<T 
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and e ^ 1/2 but close to 1/2 if 0{t) e L2{n,V{{-AY/^)), <t <T. The proof is 
similar for that in Lemma 13.21 

Remark 2.6. It is important to notice t hat f or < 7 < 1, if X p e £2(0, P((-^)^)) 
with E\\{—A)'^Xq)\\'^ < 00, Assumption 2.2 Assumption 2.3 and Assumption 
ensure the existence of the unique solution X(t) € L2(D. !?((— A)'^)) such that 



2.5 



E 



sup \\{-Arx{sm 

0<s<T 



< 00. 



In general if Xo G L2(D, X>((-A)t)) with'E\\{-A)'iXo)\\^ < 00, andO{t) e £2(0, 
then X e L2(D,2?((-^)'"™^'^'"^)) with 



E sup II 

VO<s<T 



< 00. 



More information about properties of the solution of the SPDE (1.1) can he found in 

m 



Assumption 2.7. [Initial solution Xq and extra assumption for nonlinearity 

F] 

• We assume that when initial data Xq G ^2(0, X>((-A)t)) then'E\\{-A)'i Xq\\^ <| 
00, with < 7 < 1. 

• Further we assume in our main results that for a solution X(t) of our SPDE 
(O) e L2(D,H),Vi e [0,r] with sup E||F(X(s))||2 < 00. 

0<s<T 

Remark 2.8. It is important to notice that we denote by Assumption \2.,^ i) and 
Assumption \2.^ i) to specify condition (i) of Assumption 2.3 and Assumption \2. 5| 
i £ {a, b}. 

2.1. Main results. Throughout the article we let N be the number of terms of 
truncated noise, Ijv = {1, 2, iV}'^ and take tm = mAt £ (0, T], where T = MAt for 
m,M gN. We take C to be a constant that may depend on T and other parameters 
but not on At, N or h. 

Our first result is a strong convergence result in when the non-linearity satisfies 
the Lipschitz condition of Assumption 2.3 (a). 



Theorem 2.9. Suppose that Assumptions 2.2. 2.3'^a) 2.5 (a) and (h) (with r = 
1,2^ and Assumption 2.1 are satisfied. Let X{t„i) be the mild solution of equation 
(1.1) represented by (1.3) and X^ be the numerical approximation through (2.10). 
Let 1/2 < 7 < 1 and set a — min(20,7), and let 9 S (0,1/2] be defined as in 
Assumptions^ //Aq G L2(D,X'((-^)^)) then 



(E||X(t„,) - Xtff^' < C l^t-^^'ih^ + At-) + At MAt)\ 
IfXo e L2{D,V{~A)) then 

h'' + At^^ + At\lii{At)\ + 



inf A, 



,1/2 



(E||A(U-^^lP)^'"<C 



inf A. 



-r/2N 



If we assume stronger regularity on the noise we can obtain a strong error estimate 
in the norm. 
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-1/2N 



Theorem 2.10 . Suppose that Assumptions \2.^ \2.^ h), \2.5\ a) ( with r = 2) and 
Assumption 2.7 are satisfied. Let X be the solution mild of equation (1.1) represented 
by equation (1.3). Then we have the following: 
If Xo e L2{B,V{-A)) then 

{nX{tra)-Xl^Jl,^^^Y'^<c({h + ^t^/^-H-^l^)+( inf A 

If Xn e L2iD,V{-A)) and - AXq € L2iD,H^{n)) then 

(E||X(U-X^||2, )V2<cU + AtV2-e)+/' if X 



-1/2N 



for .small enough e G (0, 1/2). 

First note in both theorems we see that if the initial data are not sufficiently smooth 
then the error is dominated by this term. This behaviour is typical of non-smooth 
error estimates. Secondly we remark that if we denote by Nh the number of vertices 
in the finite element mesh then it is well known (see for example 22 ) that if > Nh 
then 



inf A 



and 



inf A, 



-1/2 



< Ch. 



As a consequence the estimates in Theorem 2.9 and Theorem 2.10 can be expressed as 
a function of h and A< only and it is the error from the finite element approximation 
that dominates. If < then it is the error from the projection P/v of the noise 
onto a finite number of modes that dominates 

we also get an estimate in the root mean square L^(r2) norm in 



2.10 



From Theorem 

the case that the nonlinear function F satisfies Assumption |2.3| (b) . 
Finally we note that if the nodes of the finite element mesh coincide with evaluations 
of 0{x, t) then the projection operator Ph is trivial. This also leads to a computational 
advantage as we no longer need the projection, we comment further in Section [4j 

3. Proofs of main results. 

3.1. Some preparatory results. We introduce the Riesz representation oper- 
ator Rfi-.V^Vh defined by 

{-ARhV,x) = {-Av,x)=a{v,x) v&V,^x&Vh. (3.1) 

Under the regularity assumptions on the triangulation and in view of the T^— ellipticity 
(2.2), it is well known (see [7]) that the following error bounds holds 

\\RhV - v\\ + h\\RhV - v\\Hiiu) < Ch'-\\v\\„.^n), veVn r E {1, 2}(3.2) 

It follows that 

\\PhV~v\\<Ch'^\\v\\H^^n) yveVDH^in), re {1,2}. (3.3) 
We examine the deterministic linear problem. Find u € V such that 



u' ^ Au given u(0) = w, te(0,T]. 



(3.4) 
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The corresponding semi-discretization in space is : find Uh € Vh such that m'^ ~ AhUh 



where u^j-^ = PhV. The fuU discretization of (3.4) using imphcit Euler in time is given 

by 

ul+^ = (I - At PhV - Sl% Phv. 

We consider the error at t„ = rtAt and define the operator T,, from 

u{tn) -ul = {S{K) -{I -At Ah)-'' Ph)v Tr,v. (3.5) 



Lemma 3.1. The following estimates hold on the numerical approximation to (3.4-). 
Estimation in H — L^(il) norm. If v £M then 

||u(t„) - ulW = \\T,M\ < Ct-^'^ih" + At)\\v\\^ (3.6) 

and if V &V{-A) =Ur\H'^{VL) 

\\u{tr,) ~ul\\ = ||T„«|| < C{h^ + M)\\vh. (3.7) 

Estimation in H^{il) norm. If v eM then 

\\u{t^)-K\\HHn) = \mM\HHn)<C\\vUt-^/'h + t-'At). (3.8) 

Ifvev{-A) ^mnH^{n) 

h(t„) - = \\Tnv\\m(n) < CWvhih + t-^^At). (3.9) 

Finally, if v e T^{-A) and Av e II^{^) then 

|j«(<„) - = WTnvWmi^n) < C\\ - Av\\Hiin)ih + At). (3.10) 

Proof. We give here some references for the proof. Estimates in the H — L^(ri) norm 
are given in [71 [H]. In [21], A — A with Dirichlet boundary conditions. Estimates in 
the H^{n) norm are the special cases of Theorem 5.3 in [T6 where the proof is given 
for a general semi-linear parabolic problem with a locally Lipschitz nonlinear term. 
To obtain our result from [16] note that u{t) = S{t)v so that we have the analogue of 
[I6l Theorem 5.2] 

\\ut{t)\\H2in) + \\uttm<Ct-^/^v\U if vGU, 

\\utmHHn) + \\uttm<Ct-'\\v\\2 if veV{~A)^VnH\n), 

\\Mt)\\HHn) + \\uttm < Ct-^^^\\ - AvWHi^n) ^ v e V{-A) and - Ave H\n). 

Using these in the proof of [HI Theorem 5.3 ] gives the result. □ 

Our second preliminary lemma concerns the mild solution SPDE of (1.1 1. 

Lemma 3.2. Let X he the mild solution of given in (1.3), let Q < 7 < 1 and 

fl,f2G[0,T], ti<t2. 

(i) IfXo € L2{B>,V{{-A)'i)), ||(-A)("-i)/2gi/2||^^ < 00 withO< a<2 and suppose 
F satisfies Assumption \2.,'^ (a). Set a — mm{'-f, 1/2, {a — e)/2) with e e (0, a/2) very 
small enough, then 

nX{t2) - X{t,)f < C{t2 ~ t,)'^ (nXof^ + E (^^sup^ (||F(0)|| + \\X{s)\\)^ + 1 J . 
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Furthermore 

n\{X{t2) - 0{t2)) - {X{ti) - 0{ti))f < C(t2 - t^f^ < 7 < 1. 
(ii) IfXo e L2(D,2?((-A)(^+i)/2)).. ||(-A)V2gi/2||^^ < ^ 

and F{X{t)) e L2{D,H^{n)), V t e [0,r] with E ( sup ||F(X(s))||jji(o) ) < oo 

Vo<s<T / 

then 

^ 2 

2 ^ r^(-L 4- AT ( Tpll V l|2 



nnh) ' x{t,)\\' < c{t2 - t,r i euxoH^^+d + e (^^sup^iif(x(s))iui(o) ) + 1 



Proof. 

Proof of the first claim of part (i) The proof is the same as in [20j Lemma 5.10] 

with the weak assumption \\{—A)°'/^Q^''^\\hs < oo- 
Consider the difference 

X{t2)-X{ti) 

= {S{t2) - ^0 + ^ ' S{t2 - s)F{X{s))ds ~ ^ ' S{t, ~ s)F{X{s})ds 

+ (^J^ ' S{t2 ~ s)dW{s) - ' S{ti - s)dW{s) 
= 1 + 11 + 111 

so that E||X(t2) ~ ^(ii)lP < 3(E||/||2 + E||//||2 + E||mf ). 
We estimate each of the terms I, II and ///. For /, using Proposition |2.1| yields 

\\I\\ = \\S{t,){~A)-''{I--S{t2-t,)){-ArXo\\ < C{t2~t,r\\Xoh. 

Then E||/f < C{t2 - ti)^'i'E\\Xo\\^. For the term //, we have 

11= f iSit2 -s)- S{ti - s))F{X{s))ds + f ' S{t2 - s)F{X{s))ds 
Jo Jti 

= Ih+Il2. 

We now estimate each term Hi and Il2- For Hi 

i-ti 



E||//i||2 = E|| / \s{t2 -s)- S{ti - s))F{X{s))ds\\ 
Jo 



ti 



<E / \\{S{t2^3)^S{ti^s))F{Xismds 







2 



ti 



< \\iS{t2^s) 







Siti-s))\\LiLHn))d.s] b( sup \\FiX{s))\\ 

J \0<s<T 
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For < 7 < 1, Proposition |2 . 1| yields 



E||//i||' < 



\\Sih~.s)i-Ar{~A)-''il~Sit2~hm^LHn))ds) e( sup 

/ \0<s<T 

< ( r Wi-^ySih ~ - ^(^2 - ti))||L(L^(o))rfs) E { sup \\F(X{s))\\ 

< C{t2 - hf'< ( f\h^ s)-^ds) E ( sup 

\J() J \0<s<T 



\0<s<T 



<C(h-hf-'I^( sup . 

\0<s<T J 

For 7/2, using the fact that the semigroup is bounded, we have 
E||//2||2 = E|| r Sit2-s)F{X{s))dsf 

<E(^l'^\\Sih-s)FiX{smds^ 

< E f / ' < C{h - ti)'E ( sup 

\Jti J \0<s<T 



Hence, if F satisfies Assumption 2.3 (a) we have 

■TE 

\o<s<T 



E||//||2 < 2(E||//i||2 + EII//2IP) < C{t2 - ti)2^E ( sup (lli^(O)ll + \\Xis)\\) ' 



For term /// we have 



/// = f {S{t2 -s)~ S{ti - s)) dW{s) + f ' S{t2 - s)dW{s) = nil + II h- 

Jo Jti 
First using the Ito isometry property and then < a < 2 we have 

E||mi||2 = E|| /*' iSit2 - s) - Sih - s)) dW{s)\\^ 
Jo 

= / ' E|| {S{t2 ~ s) - Sih - s)) Q'^Ynsds 
Jo 

E\\{Sit2 -S)- Sih - s))(~A)-("-1)/2(-A)("-1)/2q1/2||2^^^^_ 



Using Proposition 2.1 that ||(-A)("-i)/2q1/2|| 

Hs < 00 and boundedness of S yields 

E||///i|p < r ||(-A)-("-i)/2(5(t2 -s)^ Sih - mliLHn))ds 
Jo 

\\i~A)('-^y'Sih - - Sit2 - hml^mn))ds 

rti 



<cit2-hT^' j ih-sf-'ds 
<cit2~hT-\ 
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with e G (0,a/2) small enough. 

Let us estimate E||///2||. Using the Ito isometry again, and that for < a < 2 we 
assume ||(-A)'""^^/^Q^/^||ffs < oo then WQ'^^'^Whs < oo (by taking a = 1), with 
boundedness of S yields 



= E|| 



Hence 



S{t2-.s)dW{s)\\' 



2a 



E||///||2 < 2(E||///i||2 + EII///2IP) < C(<2 - ii)2n"n(7,("-e)/2,l/2) ^ _ 

with e e (0,a/2) small enough. Combining our estimates of E||/|p, E||//|p and 
E||7//||2 ends the first part of the first claim in the lemma. 

The proof of the second claim of part (i) and part (ii) can be found in [221 Lemma 
5.10]. □ 



3.2. Proof of Theorem 



2.9 



We now estimate (E||X(i,„) - X'^\\'^)^''^ . Again 
we look at the difference between the mild solution and our numerical approximation 
( |2.10[ ). By construction of the approximation from ( |2.9| and (2.81 we have that 

X{tm) — = ^{tm) + 0{tjn) — 

= X{trn) + 0{t^) ~ [Zt + PhPNO{t^)) 

= {X{tm) ~ Z'f^) + {PN{0{t^)) - PhPN{0{t^))) + {0{tra) " Piv(0(im))) 
= 1 + 11 + 111, (3.11) 

where X{t) is given by 

X{t) = S{t)Xo + [ S{t- s)F{X{s))ds, t e [0, T]. 
Jo 



and by ^ 



Then {E\\X{t„,)~X!^,fy^^ < (E||/||2)^/' + (E||//||2) + (E||///||2)^'^ and we 
estimate each term. Since the first term will require the most work we first estimate 
the other two. 

Let us estimate (E||//|p) . To do this we use the finite element estimate (3.3 1, the 
regularity of the noise from Assumption 2.5 (a) and the fact that Pjv is bounded. 
Then for r € {1, 2} we have 



,1/2 



,1/2 



E||//|p < C/i2'-E||PAr(O(t„0)lk^(o) < Ce''E\\0{t„MH^in)- 
= \\{-AY^'^.\\ in X>((-A)''/2) and the Ito isometry yields 



Using 

= Ch^"^ 



Sit^-s)dW{s)f 
\\(^AY^'S{t^-s)Q'^'\\lsds 



-Xit„ 



ch 



2r 



(r-l)/2^1/2| 



HS- 
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Thus, since the noise is in -ff, using relation (2.12) of Assumption 2.5 yields (E||//|p)^/^ <| 
For the third term /// 



E||///||2 = E||(I - Pjv)0(t„,)|p = E||(I - Pj,){-Ar-'^{-AY'^0{t,^)\\ 



and so 



nilir < \\{l-PN){~A)-^^YL(LHn))m-Ay^'0{t„,)r < C ( .JnL A, 



We now turn our attention to the first term E||/||^. Using the definition of Sh,At in 
(2.7) we can write (2.8) as 

m— 1 
k=0 

Then using the definition of T„j from ( |3.5[ ) the first term / can be expanded 

m-l „t 



/ = T^Xo +Y. ^(^"^ - ^)^(^(^)) - ^["'At PhF{Zj: + PhPNO{tu))ds 

/c=0 

= T„,Xo +J2 S]l'Xt''^Ph{F{X{tk)) - F{Z'^ + PhPNOihWds 

k=o 

E / sl:^^,'^P,{F{X{s))-F{Xit,)))ds 



k=0 



E / iS{tr.~h)-Sl"'^;^P,)F{X{s))ds 



™-l /-ifc+i 

E / - s) - - ife))^^(X(s))ds 

7. — n 



fc=o 

= /i + /2 + h + h + Ii 



(3.12) 



Then 



For h, if Xq e L2(D,X>((-yl)T)), 1/2 < 7 < 1, equation (|3^ of Lemma [sl] gives 

mhry^'< c{t-^^'{h^+ At)) inxoiiif' 



and if € ^2(10', A)), equation (3.7) of Lemma |3.1| gives 

iMhf)'/' <Cih' + At) (EllXoll^)'/'. 

If satisfies Assumption 2.3 (a), then using the Lipschitz condition, triangle inequal- 
ity and the fact that S^^^ and Ph are an bounded operators, we have 



(EII/2IP)'/' < CE / {nF{X{tu))-F{{Zt + PHPNO{tu))f)"^ds 

k=0 

<cE / {nxitk)-x'^r)'^'ds. 
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Let us estimate (E|j/3|p) . We add in and subtract out 0{s) and 0{tk) 

/ m—l t \ ^/•^ 

(Eii/3ip)^/' = 1^"' si:::,-'^p,iFixis))-Fixit,)))ds\\A 



< (e||^ I^^^' si'X'^P,{FiXis))-FiXit,) + Ois)-Oit,mds\\A 



fe=o ■ 
"1-1 rtk+i 



+ [nT. s\:Xt^Pu {F{X{tu) + 0{s) - 0{t,)) - F{X{t,)))) dsf 

\ k=a •^*'= 

Applying the Lipschitz condition in Assumption |2.5( using the fact the semigroup is 



bounded and according to Lemma 3.2 for Xq G L2{V),'D{{—A)'^)), < 7 < 1 we 
therefore have 

imirf' < / (Eii(x(,s) - o{s)) - (x{t,) - o{t,))rf'ds 

k=o •^*'= 

< / {s-tk^ds < CAf'. 

1. — n 



1 f? 

Let us now estimate E . The analysis below follows the same steps as in 

[13] although the approximating semigroup Sh,At is different here. Applying a Taylor 
expansion to F gives 

E(||/||P)'^'</f +/f + /f, 

with 

X 1/2 





m— 1 


rtk+1 


E| 


E 








'tk 




m— 1 


ftk+i 


E| 


E 








'tk 




m— 1 


rtk+i 


E| 


E 








'tk 



V fc=0 

^1/2 

/f =(e||5: / ' / G(l - r)drd,s|p 



1/2 



where G is the expression G := P/ii^"(A:(<fc)) + r(0(s) - 0(tfe))(Cl(s) - O(tfe), 0(s) - 
O(tfc)). Using the fact that 0(^2) - 5(^2 - ii)O(ii), < ii < ta < T is independent 
of !Fti , one can show, as in [13] . that 

^E|| / si'X'^P,F'{X{t,)){0{s)~S{s-t,)0{t,))dsr. 
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Therefore as Sh.At is bounded we have 

^ (E\\si%''^PkF\X{t,)){0{s) - S{s - h}0{t,))r) ds 

V*;=0 ^-^^fc 



Using Holder's inequality, the following inequality holds 

fix)dx] <{b~a) fixfdx, (3.13) 



by assuming that / and are integrable in the bounded interval [a,b]. Using (3.131 
with Assumption 2.3 |a) and Proposition 2.1 yields 



Cm— 1 „ 
E / 
k=o -^t' 



mPhF'{X{tk))iOi.s) - S{.s - tk)Oitk))fds 

1/2 



< ^ / E||(0(s) - 5(s - t,)0{tk)Wds 



,k=0 •'^i' 



<CM^'^[Y, / ((E|10(.)~0(t,)|p)'/' + (E|l(5(.-i,)-I)0(tfe))|p)'/') ds 

\A;=0 y 

Using Assumption |2.5[ b) yields 
Let us estimate I^^. 

r*i'+^ / r , \ , , \ 1/2 



1/; 



ds 



k=0 " "^fc 



Since Phi—A)^/"^ = (— ^/i)^/^ and 5'/i_At satisfies the smoothing properties analogous 
to S{t) independently of h (see for example [161 [5]). and in particular 

we therefore have 

/f < C^E /''^' - ^fe)"'^' (E||(-^)-i/2F'(X(tfc))((5(s - h) - l)0(ife))|p) ds. 
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The usual identification oi H — L (fi) to its dual yields 

k=0 



,. n. Jtl, 



X I E f sup \{v, {-A)-^/^F'{X{tkmS{s - tk) - l)0{tu)))\\ \ ds, 



t)||<i 



where (, ) = (, ) and we change the notation merely to emphasize that H is identified to 
its dual space. The fact that (-A)''^/'^ is self-adjoint implies that {{-A)-^^^F'{X))* =| 
F'{X)*{-A)-'^/'^. This combined with the fact that V{{-Ay^'^)c H thus H = H* d 



2?((-yl)"i/2) = (X>((-A)i/2)) continuously and Assumption |2.3| yields 

I n ^ th 



X (^supKF'iX{tk)r{-A)-'/^v, {S{s - tk) - l)0{tk))\] ^ ds 
'"-1 /-tfc+i 



fe=0 •^*'= 

s 2\ 1/2 

X |e( sup ||i^'(X(ife))*(-A)-i/2t;||i||(^(s-tfc)-I)0(ifc))||_i ) I ds. 

\\v\\<l 



We also have 



/I' < - tk)-'^' (e (1 + ll^(ifc)lli)' II {S{s - tk) - l)0{tk)\\-i)') ds 

rtk+i . \ 1/4 / \l/4 

<C^ / (<„-tfe)-^/'(E(l + |lX(tfc)|li)^) (E(||5(s-tfe)-I)0(tfe)|l_i)^) ds. 
fc=o 

771 — 1 pik-\-i 1/4 

(t„,-tfc)-^/'(l+(E||X(<,)||t)'/*) / (E(|15(s-tfc)-I)0(tfe)|l_i)^) ds 

J n -'iA: 



fc = 



<CJ2 itn.^tk)-'^' / ||(-A)-('-/2+V2)(^(,_i,)_I)||^(^,(^^^^(E||0(t,.)||4)l/4ds 
fc=0 **= 

< (t™ - tfc)^'/' / ||(-A)i/2"'-/2(_A)-i(5(s - tk) - I)\\^LHn))ds. 

J n -^it- 



Using Proposition 2.1 and the fact that {—A)^^^ is bounded as r € {1,2} yields 

k^O '^^^ k^O 



771-1 rtk + 1 771-1 
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As the sum above can be bounded by 2M^/^ we have 

/|i + < C{At + Ati/2+») < C{At^'). 

Let us estimate /|'^. Usmg the fact that 5'^"Xt'^'' is bounded for any to, k with Assump- 
tion[23]and Assumption [2l] yields (with G = P^F" {X(tk) +r{0{s) - 0(ifc)))(0(s) - 
Oitk),0{s)-0{tk))) 



/f < E / ll^lf ^IU(L^(0)) / (E||G|r)^/^drd. 

I— n "'tfc "'0 

{E\\0{s) - OitkT^y^' drds 
^^E /"""((E||OW-0(t,)||t)' 



k=o 



fc=0 



Combining + and /|'^ yields the following estimate 



0<s<T 



We now estimate I4. Using equation (3.6) of Lemma 3.1 if F{X{t)) G L2(D,H) with 
sup 'E\\F{X{s))\\l < 00 (Assumption! 



2.7), we have 



m — 1 



<C(/.2 sup {E\\FiXis))\\iy'[AtY,C/! 



/ m— 1 > 

1/2 / 



0<s<T 



fc=0 



Note that we can bound AtJ2"kld tj-l by ^VT then 



(E\\l4\fY^^<C{h^ + At)( sup E||F(X(s))||2 

Vo<s<T 



1/2 



1 /9 

Finally we estimate (E||/5|p) . Using Proposition 



2.1 



we have for < ii < i2 < 7" 



(t2-tl) 



\\Sih) - S{h)\\LiLHn)) = \\i-A)Sih)i-A)-' (I - S{h - h)) himn)) < ^ 
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then 



(EII/5 



,2\l/2 



\k=0 ■"^'= 



\\Sit^~s)-S{t^-tk)\\LiLHn))ds\ sup E||F(X(s))|l- 

\0<s<T 



™-2 „tk+i 



k=0 
m-2 



S-tk 



ds 



< C At + ^ ((m - fc - l)Ai)"^ 

\ fc=0 

< C [ At + Ai"^ (m-k- 1)"M 



(s - tfe) ds 



\ fc=o / 

Noting that the sum above is bounded by ln(M) we have 

{Mhff^ < C{At + At\\n{At)\). 
For F{X{t)) e L2(D,H) with sup E\\F{X{s)\\j < 00 (Assumption 

0<s<T 

the previous estimates for the term / yields: If Xq G L2{D,'D{~A)), 



2.7), combining 



(E||/f <c[h^ + Af" + At\ ln(A<)| + J2 I - ^fe 11')'^' ^ 



m-l „t 



IfXoGL2(D,2?((-A)^)) 



imrf <c{t;;,'/'{h' + An + At\\n{At)\+J2 I '^'^ {E\\X{t,)~Xl: 



m-l „fj^^ 



k=0 



(is 



where tr — min(20,7). 
Finahy we combine ah our estimates on /, // and /// to get (E||/||^) " , (E||// 

1/2 

and (E||///||2)''" and use the discrete Gronwah lemma to complete the proof. 



,1/2 



2.10 



We now prove convergence in H^{il) and esti- 
For the proof we follow the same steps as in previous 



3.3. Proof of Theorem 

mate (E||X(i,„) - J^Jf) 1/2^ 

section for Theorem 2.9 We now estimate (3.11) in the norm. 

The estimates of the terms // and /// follow as in Section 13.21 and we find 



(e|i///|| 



1/2 



< C inf A, 



-1/2 



We concentrate instead on estimating the first term I = I1+ 12+ 13 + 14 + 15 in (3.12) 



and estimates on Ii follow immediately from Lemma 3.1 



If F satisfies Assumption 2.3 (b), then using the Lipschitz condition, the triangle 
inequality, the fact that Ph is an bounded operator and Sh.At satisfies the smoothing 
property analogous to S{t) independently of h [16J, ie 
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we have 

m-l „t. 



k=0 
'"-1 /-tfc + i 

^ / CUmFiXitk)) F{{Zji + P,,Pr,0[tu))rY'^ds 

k=o 

'"-1 /-tfc+i , 
< ^^E / C^lin^i^^) X'k\\hin)y^'ds. 

If Xq £ L2(D, 2?(— A)), again using Lipschitz condition, triangle inequality, smoothing 



property of Sh,At, but with Lemma 3.2 



gives 



fe=0 



^C'E / Cll{nF{X{s))-F{X{t,))\\)'/'ds 
™-l /-tfc+l , 

< ^^E / Cllin\x{s) ~ x{t,)\\i,^^^Y/^ds 



nM\i + ( E^sup^||F(X(s))||^i(f,) ) + 1 



2 



< C(At)(i/2-^) I^AiEC-fc j (^EllXoll^ + (^E^sup^||^^(X(.))|U.(n) ) + 1 
e e (0, 1/4) small enough. 

As in the previous theorem, we use the fact that Ai^^^^t^^^^ < 2\/T. 

k=0 _ 

For F{X{t)) e £2(0, H) with sup E\\F{X{s))\\l < 00 (Assumption 



0<s<T 



relation (3.8 1 of Lemma 3.1 we find 



™-i /-tfc+i 

inhfm^n))'^' < E / mTm-kF{x{sm%r^^^)'^'ds 

k=0 



2.71, then by 



<CAt (E*™-> + C-fcAi) ( sup E||F(X(s))|l?' 



\ fe=0 

^m-l ,-1 



0<s<T 



Note that At X^Lo ^m-fc < ln(T/Ai) to get 

/ X 1/2 'Z^^ /-^fc+i -1/9 

(eii/4|ii^i(o)) <E / {nTm-kF{x{s))\\iY^'ds 



tk 

1/2 

2 ■ 



fc=0 

< C{h + At\n{T/At)) ( sup E,\\F{X{s))\\'i 

vO<s<T 
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Finally, using the equivalency ||.||_f/i(n) = ||(-A)i/2.|| in ^((-A)^/^), we obviously 
have for < ti < ^2 < r 

\\S{t2) - Sih)U^HHn)) < C\\i-Ar/'Sih)i-A)-' (I - 5(^2 - h)) [^(^.(o)) < C^-^^^ 



so with splitting yields 

t'm.-l .tk+i 



/m-1 „tk+i \ / 

inhWl^n))'^^ <\T. - ■'^) - - tk)\\LiHHn))ds sup E||F(X(s)|l?,, 

< C \\S{tm - s) - S{tm - tm-i)\\L{m{n))ds 

'"-2 _ rtk+1 \ 

+ 51 - ^fc - / - 

fc=0 -^^^ / 

(m-2 \ 
Ail/2 + Ail/2 ^ (m-fc-l)-^/M . 
fc=0 / 

Since the sum above can be bounded by 2 we have that (E||/5|p) ^''^ < CAt^/^. 
Combining our estimates, and using that IS.t'^^^''/'^^ In (T/ At) is bounded as At — 0, 



for F{X{t)) e L2{D,M) with sup E\\F{X(s))\\l < oo (Assumption 2.7h, we have the 

0<s<T 

following. 

If e L2(D,X'(-A)) then 

/ m-l tj^^j N 

(E||/||l,.(o))^/^ < C (;^ + AtV2-ei^i/2) + ^ / t;^V2(E||x(i,) - V2rf, 

\ fc=0 / 

liXo e L2(D,2?(-^)) with -AXq G ^2(0, iJ^l^^)) then 

/ m-l t^^i \ 

(E||/|1^.(,,))V2 < C (/i + AiV2-) + ^ / t-V2(E||x(t,) - X^\\],,^^^y/^ds . 

\ k=o -^'fc / 

where C > depending of the T, the initial solution Xq, the mild solution X, the 
nonlinear function F. 

Combining our estimates (E||/||2) (E||//|j2) and (E||///||2) and using the 
discrete Gronwall lemma concludes the proof. 

4. Numerical Simulations. We consider for our numerical simulations two ex- 
amples. Below we consider the error estimates in L'^{fl) since this is standard in many 
of the applications for porous flow and the discrete L'^{Q,) norm is straightforward to 
implement for different types of boundary conditions. 
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4.1. A linear reaction— diffusion equation. As a simple example consider the 
reaction diffusion equation 

dX = [DAX - 0.5X)dt + dW given A(0) Aq, 

in the time interval [0, T] with diffusion coefficient D = 1 and homogeneous Neu- 
mann boundary conditions on the domain — [0,Li] x [0,^2]. The eigenfunctions 
{ep^e^^''}i,j>o of the operator A here are given by 

= ef)(x) = y|cos(Af).), A(')=0, Af ^ = ^ 

where / S {1,2} , x G and i — 1,2,3, ••• with the corresponding eigenvalues 

{Ai,j}ij>o given by Xij = (Af + (A^)^. We take Li = L2 = 1. Notice that 
A = DA does not satisfy Assumption |2.2| as is an eigenvalue. During the simu- 
lations we need to manage the singularity of (4.2) at Ag = or use the perturbed 
operator A = DA -I- el, e > 0. 

Our function F{u) = — 0.5u is linear and obviously satisfies Assumption 2.3 (a) with 
V = L'^{n). In our simulation we take Xq = 0, then X{t) e L2(D,H) t e [0,r] 
(according to Remark [2^ and we obviously have F{X{t)) € L2(D,HI) t e [0,r], so 



that Assumption |2.7| is satisfied. In general we are interested in nonlinear F however 
for this linear system we can find an exact solution to compare our numerics to. 
Recall that we assumed for convenience that the eigenfunctions of covariance operator 
Q and A are the same. In order to fulful Assumption 2.5 with V — L^{n), then 6* = 1/2 



in Assumption |2.5| [b), we take in the representation (1.2) 

*j = '"^^"^^ r > 0, and 5 > small enough. (4.1) 



We obviously have 

-2 /■2 , ■2\-(l+'5) 



AU'9..<^' E {^+fy'^'<°° re {1,2} 



thus Assumption |2 . 5[ a) is satisfied. In the implementation of our modified scheme at 
every time step, 0{tk+i) is generated using 0(tk) from the following relation 



where 0(0) = 0. We expand in Fourier space and apply the Ito isometry in each mode 
and project onto N modes to obtain for fc 1, 2, . . . , M — 1 

(e„ 0(tfc+i)) - e-^'^*(e„ 0(t,)) + (l - e-^^'^*)) ' i?,,^, (4.2) 

where Rij. are independent, standard normally distributed random variables with 
means and variance 1, and i S In = {1, 2, 3, A}^. These are the linear functionals 
used by Jentzen and Kloeden in [12]. The noise is then projected onto the finite 
element space by Ph. We evaluated 0{t) at the mesh vertices, thus the projection Ph 
becomes trivial. 

In our simulations we examine both the finite element and the finite volume discretiza- 
tion in space. For the cell center finite volume discretization we take Ax — Ay = 
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1/100. The finite element triangulation was constructed so that the center of the 
control volume for the finite volume method was a vertex in finite element mesh. In 
all our simulations in this paper, we denote by 'Modifiedlmplicitfemr', r € {1,2} for 
the graphs with the finite element in space discretization and our modified scheme 
for time discretization. For the finite volume discretization we have implemented 
both the new modified method denoted by 'Modifiedlmplicitfvmr', r £ {1,2} and 
a standard implicit Euler-Maruyama method denoted Tmplicitfvmr', r G {1,2}. In 
Figure |4.1[ a), we see that the observed rate of convergence for the finite element 
discretization agrees with Theorem 2.9 with V — L'^{n) and ( 



1/2, and rate of 

convergence in At is very close to 1 for r e {1, 2} in our modified scheme. We also 
observe that the finite element and finite volume methods for space discretization give 
the same errors as the mesh is regular. More importantly we see that the error using 
the new modified scheme is smaller in all the cases than using the standard implicit 
Euler-Maruyama scheme. Indeed we observe numerically a slower rate of convergence 
for the the standard scheme of 0.65 for r = 1 and 0.98 for r = 2 respectively. This 
is compared to 0.9960 for r = 1 and 1.0074 for r = 2 with the modified scheme. We 
also observe that the error decreases as the regularity increases from r = 1 to r = 2. 




■ Implicittvml 

Implicittvm2 
~* — Modifiedlmplicitfemi 

Modifiedlmplicitfem2 ; 
— * — Modifiedlmplicitfvml 

Modifiedlmplicilfvm2 



Fig. 4.1. Convergence in the root mean square norm at T = 1 as a function of At. We 
show convergence for noise where the stochastic process O is in H'^{Q) respect to space variable, 
r S {1,2} and 5 = 0.05 in relation l |4.1[ l for finite element and finite volume space discretizations. 
We also show convergence of the standard semi-implicit scheme for the finite volume discretization. 
We used here 30 realizations. 



4.2. Stochastic advection diffusion reaction. As a more challenging exam- 
ple we consider the stochastic advection diffusion reaction SPDE ( 1.6 ), with D = 10~^ 
and mixed Neumann-Dirichlet boundary conditions. The Dirichlet boundary condi- 
tion is X ^ 1 at T ^ {{x,y) : x = 0} and we use the homogeneous Neumann 
boundary conditions elsewhere. Our goal here is to show that with the well known 
eigenvalues and eigenfunctions of the operator A with Neumann (or Dirichlet) bound- 
ary conditions, we can apply the new scheme to mixed boundary conditions for the 
operator A ~ DA without explicitly having eigenvalues and eigenfunctions of A. In- 

(see 



deed in term of the abstract setting in (1.1), using the trace operator 71 
[T7]) in Green's theorem yields 



dX = {AX + Fi{X) + b{X))dt + dW, 



(4.3) 
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where 



{Au,v) — — D\Iu\7vdx, {bu,v) ~ / jiujovda, jqv = v v e H^{ri)^4:A) 
Jn Jr 

(4.5) 



dx 



u e {x e H^(n) : ^ ^0 in TA, Ti ^ dn\r. 
ov 



In the abstract setting in ( 1.1 ), the hnear operator is yl = DA using only homogeneous 



Neumann boundary. The explicit expression of 6 is unknown. The finite volume 
method (finite element or finite difference methods) uses a natural approximation of 
b (see [ini [17] for finite volume approximation). The nonlinear term is F = Fi + 6 
where 



Fi(w) = -V-(qu)- 



{u + iy 



u e 



(4.6) 



As h is linear, F clearly satisfies Assumption 2.3 (b). We use a heterogeneous medium 



with three parallel high permeability streaks, 100 times higher compared to the other 
part of the medium. This could represent for example a highly idealized fracture 



pattern. We obtain the Darcy velocity field q by solving the system (1.7) 



with Dirichlet boundary conditions r]^ = {0, l}x [0, 1] and Neumann boundary T\i =| 
(0, 1) X {0, 1} such that 



m 
in 



{0} X [0, 1] 
{ill X [0,1] 



-fc Vp(x,t) • n = 



i AT- 



To deal with high Peclet flows we discretize in space using finite volumes. We can 
write the semi-discrete finite volume method as 



(4.7) 



where here Ah is the space discretization of DA using only homogeneous Neumann 
boundary conditions and Phh{X^) comes from the approximation of diffusion flux on 
the Dirichlet boundary condition side (see [IHlUn])- Thus we can form the noise as in 
Secti on|4.1| with the eigenvalues function of A with full Neumann boundary conditions 
and \A.l\. 



Figure 4.2(a) shows the convergence of the modified method and standard implicit 
with 0{t) € L2(D,i?''), r e {1,2}, < t < T. We observe that the temporal 
convergence order is close to 1/4 for all the schemes. The predicted order 0.5 for r = 2 



four our modified scheme is probably not achieved due to the fact that Assumption |2.7| 
is not satisfied. We do note, however, that even if the conditions of the theorem are 
not fully satisfied we still obtain an improvement in the accuracy over the standard 



semi-implicit method. Figure 4.2(b) shows the streamline of the velocity field, Figure 
4.2(c) shows the mean of 1000 realizations of the "true solution" (with the smallest 
time step At = 1/7680) for r = 1 while Figure 4.2(d) shows a sample of the "true 
solution" with r = 1. 
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Fig. 4.2. (a) Convergence of the root mean square norm at T = 1 as a function of At with 
1000 realizations with Ax = Ay = 1/100, Xg = 0. The noise is white in time and the stochastic 
process 0{t) £ H^{fl) respect to space variable, r £ {1)2} and S = 0.05 in relation ( |4.1| . The 
temporal order of convergence in time is 1/4, (c) is the streamline of the velocity field. In (b) 
we plot a sample the 'true solution' for r = 1 with At = 1/7680 while (d) is the mean of 1000 
realizations. 
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